OpenFoam Assignment | The Mach 3 Step
Abstract
–SAMPLE ABSTRACT–
Keywords: computational fluid dynamics, CFD, incompressible, Paraview, R, Python, coe347, spring 2022.
Motivation
The flow over a Mach 3 step is a well-studied problem in CFD (computational fluid dynamics), making it a good choice for convergence and performance studies. A simple set of studies are undertaken here.
Implementation
We implement all simulation with OpenFoam, analysis with Paraview and Python3, and documentation code in R Xie, Dervieux, and Riederer (2020).
Theory : Governing Equations
Mesh Assembly : Preliminaries
We assemble a preliminary solution at \(T = 4s\) with 6300 cells, for the following inlet conditions: \(M=3\),\(p=1Pa\),\(T=1K\),\(\rho = \gamma \frac{kg}{m^3}\) with \(\gamma = 1.4\). All results have been non-dimensionalized: velocity, pressure, density and temperature are all measured with respect to inlet conditions.
Mesh
Density
Velocity X,Y, Magnitude
Pressure
Mach Number
Speed of Sound
Temperature
Convergence
Three meshes were used, with increasing numbers of cells and smaller timesteps, as follows: (note the timesteps are approximate as the solver, RhoCentralFoam, adjusts timesteps as necessary to maintain a maximum Courant number of 0.2).
Run_10_64 : 6300 cells, dx = 0.04000, dt = 0.00110
Run_20_64 : 25200 cells, dx = 0.02000, dt = 0.00056
Run_32_64 : 64512 cells, dx = 0.00625, dt = 0.00034
Improvements in Pressure, Density, and Y-Velocity with increasing resolution:
Note the upper shock turbulent features (the small waves), and how they resolve much more appropriately with higher resolution meshes.
Pressure (from Run_10 to Run_32):
Density:
Y-velocity:
Wall Distributions w.r.t Resolution Levels:
To examine the convergence of solution fields with a refinement of the mesh and associated reduction in grid spacing, plots of the pressure, temperature, fluid density, and x and y components of velocity were developed against the distance along four surfaces of interest, namely, the bottom wall, “A”, the vertical step wall, “B”, the bottom step wall, “C”, and the top wall, “D”. The purpose of generating these plots was to observe the underlying patterns, and convergence of these solution fields to theoretical behavior.
Wall A Pressure, Density, Temperature:
In referencing figures 1-3 for the pressure, density, and temperature fields across the bottom wall, “A”, it is shown that the profiles are very similar, namely, the properties are approximately constant from the beginning of the wall to about 0.3 meters where they experience a sharp gradient as the flow moves across the normal shock wave boundary where theoretical behavior suggests that fluid density, pressure, and temperature increase rapidly across shock waves. These sharp gradients effectively confirm this expected behavior. Shortly after these sharp gradients (just removed from 0.3 meters from the beginning of the wall), the properties level off and nearly reach steady state. This behavior agrees with theory as the region downstream of the flow is isentropic (constant entropy) , and in this region, all thermodynamic properties are expected to remain relatively constant, or increase marginally with position. Furthermore, from all three plots we note that the solutions for the different spatial resolutions nearly overlap over [0,0.3 meters], and over [0.32,0.6 meters], while the largest differences in values are apparent over the normal shock wave region as expected. Additionally, in referencing the pressure plot, it is shown that the stagnation pressure approaches the theoretical value of approximately 12.06 Pa near the end of the bottom wall, where the simulation with the greatest refinement (h = 0.00625) most closely approximates this expected value.
Wall B Pressure, Density, Temperature:
In referencing figures 4-6, the pressure, density, and temperature solutions all exhibit very similar behavior with position across the vertical step wall. In particular, the profiles are approximately inverted parabolic relationships where the properties decrease from a maximum value at the stagnation point to a minimum value at the top of the step wall where an expansion fan is realized.
Wall C Pressure, Density, Temperature:
In referencing figures 7-9, it is again shown that the three solution fields exhibit similar properties with respect to position across the horizontal step wall. As expected, there is a large gradient at approximately 0.7 meters from the left end of the wall where the flow exhibits shock reflection at the intersection of two oblique shock waves.
Wall D Pressure, Density, Temperature:
Shifting gears to figures 10-12, it is observed that the solution fields again remain relatively constant over 0.6 meters before increasing sharply where the flow exhibits its first shock reflection point. It then levels off between 0.7 meters and roughly 2.25 – 2.4 meters before experiencing a second sharp gradient as the flow encounters a second shock reflection point, before levelling off again. These behaviors thus agree with theoretical predictions.
Velocity X and Y components along A,B,C,D :
Transitioning now to the velocity field solutions, it is shown in figure 13 that the x-component of velocity remains constant before dropping off sharply at the shock wave where the Mach number decreases, and ultimately decreases linearly until the stagnation point. Figure 14 shows that the y-component of the velocity also begins at a constant value until reaching the shock wave where the velocity increases somewhat sharply before exhibiting a linear relationship. Interestingly, the solution with the largest spatial resolution exhibits the steepest gradient from the shockwave onwards, while the one with the smallest resolution approximately levels out.
Looking at figure 15, it is shown that the x-component of velocity begins at a constant value before increasing roughly linearly. Interestingly, the spatial resolution greatly impacts the location where the velocity begins to increase. Conversely, figure 16 shows that the y-component of velocity increases roughly parabolically with position, and there is nearly an overlap in the solutions for each of the resolutions. These behaviors are reasonable given that as position along the vertical step wall is increased, the flow moves further from the shock wave region (where the velocity is small given low Mach Number).
Looking now at figure 17, it is shown that the x-component of velocity increases rapidly at the beginning of the wall near the Mach Reflection region before levelling off, and then dropping sharply at the horizontal position corresponding to the first shock reflection point. It then increases from the first to the second shock reflection point. In referencing figure 18, it is shown that the y-component of velocity drops off significantly at the beginning of the wall before reaching negative values in what is likely a small recirculation region. It then levels off at approximately 0 around the first shock reflection point.
In transitioning to figure 19, it is observed that the x-component of velocity drops off sharply at the point at the shockwave along the upper wall before gradually increasing to a local maximum at the second shock reflection point, and then dropping off sharply. Lastly, in looking at figure 20, it is observed that the y-component of velocity is approximately zero over the entire upper wall until reaching the second shock reflection point.
Stagnation Pressure at Step
Note from the first two pressure plots, the pressure at the base of the step is slightly over 12Pa. The simulation measures a value of about 12.05Pa, which is in agreement with the theoretical calculation performed earlier.
Parallel Performance Improvements
Define the speedup as the time taken for one processor (\(T(1)\)) divided by the time taken for \(P\) processors (\(T(P)\)). Similarly, let the parallel efficiency \(\eta\) be the speedup divided by \(P\).
Note that the tables report the total \(T(P)\), not that of single steps. The single step \(T(P)\) is used instead for the plots.
On a mesh of 64512 cells, for \(P\) ranging from 1 to 64 in steps of 2, the following scaling is obtained.
\(T(P)\) Speedup \(\eta\)
4870.28 1.0000 1.0000
2745.57 1.7738 0.8869
1493.30 3.2614 0.8154
0865.21 5.6290 0.7036
0513.54 9.4837 0.5927
0337.99 14.409 0.4503
0247.94 19.642 0.3069
A slightly improved scaling is obtained for a mesh of 1000188 cells.
\(T(P)\) Speedup \(\eta\)
270763.66 1.0000 1.0000
034192.10 2.0239 1.0119
034129.40 3.9783 0.9946
039850.40 7.6188 0.9524
021542.00 14.094 0.8808
012256.70 24.765 0.7739
008170.29 37.171 0.5808
Clearly the larger mesh case scales much more nicely, which can be expected due to the more even load distribution.
References
Appendix
Thank you so much for reading this work!